Analysis date: 2022-08-07

Setup

Load libraries

library(plyr)
library(gtools)
library(pheatmap)
library(Matrix)
library(Hmisc)
library(ggpubr)
library(ggbeeswarm)
library(DESeq2)
library(tidyverse)
library(vsn)
library(fdrtool)
library(limma)
library(Rtsne)
library(MultiAssayExperiment)
library(gplots)
library(viridis)

Load data

source("data/Figure_layouts.R")
## Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
## `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
load("data/CLL_Proteomics_Setup.RData")
load("data/CLL_Proteomics_ConsensusClustering.RData")

load("data/GOterm_BCR.RData")
load("data/GOterm_Spliceosome.RData")
new_colData <- left_join(colData(multiomics_MAE) %>% as_tibble(),
  enframe(CCP_group5, value = "PG", name = "patient_ID") ) 
## Joining, by = "patient_ID"
new_colData <- left_join(new_colData,
  enframe(CCP_group6_RNA, value = "CCP6_RNA", name = "patient_ID") ) %>%
  as.data.frame() %>%
  mutate(rn = patient_ID) %>%
  column_to_rownames(., var = "rn")
## Joining, by = "patient_ID"
colData(multiomics_MAE)$PG <- as.factor(new_colData$PG)
colData(multiomics_MAE)$CCP6_RNA <- as.factor(new_colData$CCP6_RNA)

Analysis

Dimension reduction plots

tsne (Barnes-Hut-SNE)

t-SNE colored by genetic alterations and groups

set.seed(10)
rtsne_out <- Rtsne( t(assay(multiomics_MAE[prot_few_nas , ,"proteomics"])), perplexity = 10 )
## Warning: 'experiments' dropped; see 'metadata'
rtsne_out_df <- rtsne_out$Y %>% as.data.frame() %>% as_tibble()
rtsne_out_df$pat_ID <- colnames(assay(multiomics_MAE[prot_few_nas , ,"proteomics"]))
## Warning: 'experiments' dropped; see 'metadata'
rtsne_out_df <- left_join(rtsne_out_df, 
                          wideFormat(multiomics_MAE[, ,c("SNPs","chrom_abber", "health_record_bin") ], colDataCols = c("gender", "treatment_status", "doehner_groups", "PG")) %>% as_tibble(), 
                          by=c("pat_ID"="primary") )
## Warning: 'experiments' dropped; see 'metadata'
## Warning: 'ExperimentList' contains 'data.frame' or 'DataFrame',
##   potential for errors with mixed data types
rtsne_out_df <- mutate_at(rtsne_out_df, colnames(rtsne_out_df %>% dplyr::select(SNPs_ATM:health_record_bin_treated)), as.logical)
#rtsne_out_df <- rtsne_out_df %>% replace(is.na(.), "unknown")

message("t-SNE colored by consensus cluster groups")
tsne_CCP_P_plot <-  ggplot(rtsne_out_df, aes(V1, V2)) +
  geom_point(aes(fill=PG), color="grey", shape=21) +
  scale_fill_manual(values = colors_CCP)+
  pp_sra +
  guides(color=guide_legend(title="PG")) 
tsne_CCP_P_plot+  theme(aspect.ratio=1, plot.title = element_text(size = 30))

message("t-SNE colored by combination of IGHV status and trisomy12")
tsne_IGHVtris_plot <- ggplot(rtsne_out_df %>% filter(health_record_bin_IGHV_mutated!="unknown", chrom_abber_trisomy12!="unknown"), aes(V1, V2)) +
  geom_point(aes(fill=interaction(health_record_bin_IGHV_mutated, chrom_abber_trisomy12)), color="grey", shape=21) +
  scale_fill_manual(values = colors_CCP[c(1,3,4,2)] )+
  pp_sra_noguides 
tsne_IGHVtris_plot +  theme(aspect.ratio=1, plot.title = element_text(size = 30)) #+

sapply(colnames(rtsne_out_df)[c(4,5,8:22)], function(var){
  print(var)
  print(ggplot(rtsne_out_df, aes(V1, V2)) +
          geom_point(aes_string(color=var), size=3) +
    scale_color_manual(values = c( "#0571b0", "#ca0020", "grey"))  + 
          ggtitle(gsub("SNPs_", "", gsub("chrom_abber_", "", gsub( "health_record_bin_", "", var) ))) +
      pp_sra_noguides +
      theme(aspect.ratio=1, plot.title = element_text(size = 30)) )
})
## [1] "gender"

## [1] "treatment_status"

## [1] "SNPs_ATM"

## [1] "SNPs_BIRC3"

## [1] "SNPs_EGR2"

## [1] "SNPs_NOTCH1"

## [1] "SNPs_POT1"

## [1] "SNPs_SF3B1"

## [1] "SNPs_TP53"

## [1] "SNPs_XPO1"

## [1] "chrom_abber_del11q"

## [1] "chrom_abber_trisomy12"

## [1] "chrom_abber_del13q14"

## [1] "chrom_abber_del17p13"

## [1] "chrom_abber_gain8q24"

## [1] "health_record_bin_komplex_abberant_karyotype"

## [1] "health_record_bin_treated"

##             gender           treatment_status SNPs_ATM         SNPs_BIRC3      
## data        tbl_df,26        tbl_df,26        tbl_df,26        tbl_df,26       
## layers      list,1           list,1           list,1           list,1          
## scales      ScalesList,2     ScalesList,2     ScalesList,2     ScalesList,2    
## mapping     uneval,2         uneval,2         uneval,2         uneval,2        
## theme       theme,93         theme,93         theme,93         theme,93        
## coordinates CoordCartesian,5 CoordCartesian,5 CoordCartesian,5 CoordCartesian,5
## facet       FacetNull,2      FacetNull,2      FacetNull,2      FacetNull,2     
## plot_env    ?                ?                ?                ?               
## labels      list,4           list,4           list,4           list,4          
## guides      list,4           list,4           list,4           list,4          
##             SNPs_EGR2        SNPs_NOTCH1      SNPs_POT1        SNPs_SF3B1      
## data        tbl_df,26        tbl_df,26        tbl_df,26        tbl_df,26       
## layers      list,1           list,1           list,1           list,1          
## scales      ScalesList,2     ScalesList,2     ScalesList,2     ScalesList,2    
## mapping     uneval,2         uneval,2         uneval,2         uneval,2        
## theme       theme,93         theme,93         theme,93         theme,93        
## coordinates CoordCartesian,5 CoordCartesian,5 CoordCartesian,5 CoordCartesian,5
## facet       FacetNull,2      FacetNull,2      FacetNull,2      FacetNull,2     
## plot_env    ?                ?                ?                ?               
## labels      list,4           list,4           list,4           list,4          
## guides      list,4           list,4           list,4           list,4          
##             SNPs_TP53        SNPs_XPO1        chrom_abber_del11q
## data        tbl_df,26        tbl_df,26        tbl_df,26         
## layers      list,1           list,1           list,1            
## scales      ScalesList,2     ScalesList,2     ScalesList,2      
## mapping     uneval,2         uneval,2         uneval,2          
## theme       theme,93         theme,93         theme,93          
## coordinates CoordCartesian,5 CoordCartesian,5 CoordCartesian,5  
## facet       FacetNull,2      FacetNull,2      FacetNull,2       
## plot_env    ?                ?                ?                 
## labels      list,4           list,4           list,4            
## guides      list,4           list,4           list,4            
##             chrom_abber_trisomy12 chrom_abber_del13q14 chrom_abber_del17p13
## data        tbl_df,26             tbl_df,26            tbl_df,26           
## layers      list,1                list,1               list,1              
## scales      ScalesList,2          ScalesList,2         ScalesList,2        
## mapping     uneval,2              uneval,2             uneval,2            
## theme       theme,93              theme,93             theme,93            
## coordinates CoordCartesian,5      CoordCartesian,5     CoordCartesian,5    
## facet       FacetNull,2           FacetNull,2          FacetNull,2         
## plot_env    ?                     ?                    ?                   
## labels      list,4                list,4               list,4              
## guides      list,4                list,4               list,4              
##             chrom_abber_gain8q24 health_record_bin_komplex_abberant_karyotype
## data        tbl_df,26            tbl_df,26                                   
## layers      list,1               list,1                                      
## scales      ScalesList,2         ScalesList,2                                
## mapping     uneval,2             uneval,2                                    
## theme       theme,93             theme,93                                    
## coordinates CoordCartesian,5     CoordCartesian,5                            
## facet       FacetNull,2          FacetNull,2                                 
## plot_env    ?                    ?                                           
## labels      list,4               list,4                                      
## guides      list,4               list,4                                      
##             health_record_bin_treated
## data        tbl_df,26                
## layers      list,1                   
## scales      ScalesList,2             
## mapping     uneval,2                 
## theme       theme,93                 
## coordinates CoordCartesian,5         
## facet       FacetNull,2              
## plot_env    ?                        
## labels      list,4                   
## guides      list,4
message("There is one trisomy12 negative patient which clusters with all of the other trisomy12 patients. Does it have a subclonal mutations?")
## There is one trisomy12 negative patient which clusters with all of the other trisomy12 patients. Does it have a subclonal mutations?
metadata(multiomics_MAE)$fish_df_clonsizes["trisomy12" ,] %>% 
  dplyr::select(rtsne_out_df %>% filter(chrom_abber_trisomy12==FALSE) %>% arrange(desc(V1)) %>% slice(1:2) %>% .$pat_ID)

BCR signaling

Boxplots mean abundance

BCR_genes_mean <- assay(multiomics_MAE[BCR_genes, ,"proteomics"]) %>% colMeans(na.rm = TRUE) %>% enframe()
## Warning: 'experiments' dropped; see 'metadata'
## harmonizing input:
##   removing 522 sampleMap rows not in names(experiments)
##   removing 13 colData rownames not in sampleMap 'primary'
left_join(BCR_genes_mean, 
          wideFormat(multiomics_MAE[ c("IGHV_mutated", "trisomy12"), ,]) %>% as_tibble(), 
          by=c("name"="primary")) %>%
  mutate_at(.vars =c("chrom_abber_trisomy12", "health_record_bin_IGHV_mutated" ), as.logical ) %>%
  filter(!is.na(chrom_abber_trisomy12)) %>%
  ggplot(aes( chrom_abber_trisomy12, value ))+
  geom_boxplot(aes(fill=chrom_abber_trisomy12)) + geom_beeswarm() + 
  stat_compare_means(method = "t.test", label = "p.signif", comparisons = list(c("TRUE", "FALSE"))) + 
  ggtitle("BCR protein abundance") + 
  scale_fill_manual(values=c("#0571b0", "#ca0020")) +
  pp_sra_nox +
  ylab("Mean abundance of BCR proteins")+
  theme(aspect.ratio=2) +
  theme(legend.position = "bottom")
## Warning: 'ExperimentList' contains 'data.frame' or 'DataFrame',
##   potential for errors with mixed data types
## Warning: 'ExperimentList' contains 'data.frame' or 'DataFrame',
##   potential for errors with mixed data types

## Warning: 'ExperimentList' contains 'data.frame' or 'DataFrame',
##   potential for errors with mixed data types
## harmonizing input:
##   removing 428 sampleMap rows not in names(experiments)

left_join(BCR_genes_mean, 
          wideFormat(multiomics_MAE[ c("IGHV_mutated", "trisomy12"), ,]) %>% as_tibble(), 
          by=c("name"="primary")) %>%
  mutate_at(.vars =c("chrom_abber_trisomy12", "health_record_bin_IGHV_mutated" ), as.logical ) %>%
  filter(!is.na(health_record_bin_IGHV_mutated)) %>%
  ggplot(aes( health_record_bin_IGHV_mutated, value ))+
  geom_boxplot(aes(fill=health_record_bin_IGHV_mutated)) + geom_beeswarm() + 
  stat_compare_means(method = "t.test", label = "p.signif", comparisons = list(c("TRUE", "FALSE"))) + 
  ggtitle("BCR protein abundance") + 
  scale_fill_manual(values=c("#0571b0", "#ca0020")) +
  pp_sra_nox +
  ylab("Mean abundance of BCR proteins")+
  theme(aspect.ratio=2) +
  theme(legend.position = "bottom")
## Warning: 'ExperimentList' contains 'data.frame' or 'DataFrame',
##   potential for errors with mixed data types

## Warning: 'ExperimentList' contains 'data.frame' or 'DataFrame',
##   potential for errors with mixed data types

## Warning: 'ExperimentList' contains 'data.frame' or 'DataFrame',
##   potential for errors with mixed data types
## harmonizing input:
##   removing 428 sampleMap rows not in names(experiments)

BCR_CCP_P_plot <- left_join(BCR_genes_mean, 
          colData(multiomics_MAE) %>% as_tibble() %>% select(name=patient_ID, PG ),
          by=c("name")) %>%
  filter(!is.na(PG)) %>%
  ggplot(aes( PG, value ))+
  geom_boxplot(aes(fill=PG)) + geom_beeswarm() + 
  ggtitle("BCR protein abundance") + 
  scale_fill_manual(values=colors_CCP) +
  pp_sra +
  ylab("Mean abundance of BCR proteins")+
  #stat_compare_means(label = "p.signif", method = "t.test",
  #                   ref.group = ".all.", label.y = 0.2, hide.ns = TRUE) +
  guides(fill=guide_legend(title="PG"))
BCR_CCP_P_plot + theme(aspect.ratio=1) + 
  stat_compare_means(method = "anova", label.y = 0.22, hjust=0)

Heatmaps

Selected BCR proteins

sel_BCR <- c("ZAP70",  "IGHM", "CD79A", "CD79B",  "SYK", "PLCG2",  "BTK",  "PTPN6", 
             "MAPK1",  "PIK3CD", "AKT1", "IKBKB")

tmp_BCR <- wideFormat(multiomics_MAE[sel_BCR, ,"proteomics"], colDataCols = c("PG", "IGHV_mutated" )) %>% 
  as_tibble()
## Warning: 'experiments' dropped; see 'metadata'
## harmonizing input:
##   removing 522 sampleMap rows not in names(experiments)
##   removing 13 colData rownames not in sampleMap 'primary'
tmp_BCR_mx <- tmp_BCR %>% dplyr::select(-primary, -PG, -IGHV_mutated) %>% as.matrix()
rownames(tmp_BCR_mx) <- tmp_BCR$primary
colnames(tmp_BCR_mx) <- gsub("proteomics_", "", colnames(tmp_BCR_mx))
tmp_BCR_anno <- tmp_BCR[, c("PG", "IGHV_mutated")] %>% as.data.frame()
rownames(tmp_BCR_anno) <- tmp_BCR$primary
tmp_BCR_anno$IGHV[tmp_BCR_anno$IGHV_mutated==0] <- "U-CLL"
tmp_BCR_anno$IGHV[tmp_BCR_anno$IGHV_mutated==1] <- "M-CLL"
tmp_BCR_anno <- tmp_BCR_anno %>% select(-IGHV_mutated)

ann_colors = list(
  PG=c("1"= colors_CCP[1], "2"= colors_CCP[2], "3"= colors_CCP[3], "4"= colors_CCP[4], "5"= colors_CCP[5], "6"= colors_CCP[6] ),
  IGHV=c("U-CLL"= "#0571b0", "M-CLL"= "#ca0020"))

tmp_BCR_anno$PG <- factor(tmp_BCR_anno$PG, levels = c(5,1:4,6))

breaks= seq(min(tmp_BCR_mx), max(tmp_BCR_mx), 0.1)^2
breaks= sort(c(-breaks, breaks))
breaks <- breaks[! (breaks < min(tmp_BCR_mx) | breaks > max(tmp_BCR_mx) )]

pat_order_hclust <- sapply(c(5,1:4,6), function(P){
  hc <- hclust(dist(tmp_BCR_mx[rownames(tmp_BCR_anno[tmp_BCR_anno$PG==P,]), ] ))
  hc$labels[hc$order]
}) %>% unlist

PG_BCR_proteins_pheatmap <- pheatmap(t(tmp_BCR_mx[pat_order_hclust, ]), 
         annotation_col = tmp_BCR_anno, annotation_colors = ann_colors, scale = "row", cluster_cols = FALSE, 
         color = inferno(length(breaks)), border_color = NA, 
         gaps_col = (which(!tmp_BCR_anno %>% rownames_to_column() %>% arrange(PG, IGHV) %>% .$PG %>% duplicated())-1)[-1],  
         breaks = breaks , cutree_rows = 4, show_colnames = F, treeheight_row = 0, fontsize_row = 8)

PG_BCR_proteins_pheatmap_notscaled <- pheatmap(t(tmp_BCR_mx[pat_order_hclust, ]), 
         annotation_col = tmp_BCR_anno, annotation_colors = ann_colors, cluster_cols = FALSE, 
         color = inferno(length(breaks)), border_color = NA, 
         gaps_col = (which(!tmp_BCR_anno %>% rownames_to_column() %>% arrange(PG, IGHV) %>% .$PG %>% duplicated())-1)[-1],  
         cutree_rows = 4, show_colnames = F, treeheight_row = 0, fontsize_row = 8)

All KEGG BCR signaling proteins

tmp_BCR <- wideFormat(multiomics_MAE[BCR_genes, ,"proteomics"], colDataCols = c("PG", "IGHV_mutated" )) %>% 
  as_tibble()
## Warning: 'experiments' dropped; see 'metadata'
## harmonizing input:
##   removing 522 sampleMap rows not in names(experiments)
##   removing 13 colData rownames not in sampleMap 'primary'
tmp_BCR_mx <- tmp_BCR %>% dplyr::select(-primary, -PG, -IGHV_mutated) %>% as.matrix()
rownames(tmp_BCR_mx) <- tmp_BCR$primary
colnames(tmp_BCR_mx) <- gsub("proteomics_", "", colnames(tmp_BCR_mx))
tmp_BCR_anno <- tmp_BCR[, c("PG", "IGHV_mutated")] %>% as.data.frame()
tmp_BCR_anno$PG <- factor(tmp_BCR_anno$PG,levels = c(5,1:4,6)) #
rownames(tmp_BCR_anno) <- tmp_BCR$primary

ann_colors = list(
  PG=c("1"= colors_CCP[1], "2"= colors_CCP[2], "3"= colors_CCP[3], "4"= colors_CCP[4], "5"= colors_CCP[5], "6"= colors_CCP[6] ))

message("Heatmap of abundance BCR signaling proteins ordered according to PG")
## Heatmap of abundance BCR signaling proteins ordered according to PG
pheatmap(tmp_BCR_mx[(tmp_BCR_anno %>% rownames_to_column() %>% arrange(PG, IGHV_mutated) %>% .$rowname), ], 
         annotation_row = tmp_BCR_anno, annotation_colors = ann_colors, scale = "column", cluster_rows = FALSE, color = inferno(10), border_color = NA, 
         gaps_row = (which(!tmp_BCR_anno %>% rownames_to_column() %>% arrange(PG, IGHV_mutated) %>% .$PG %>% duplicated())-1)[-1])

Spliceosome

All KEGG proteins

splice_genes_mean <- assay(multiomics_MAE[splice_genes, ,"proteomics"]) %>% colMeans(na.rm = TRUE) %>% enframe()
## Warning: 'experiments' dropped; see 'metadata'
## harmonizing input:
##   removing 522 sampleMap rows not in names(experiments)
##   removing 13 colData rownames not in sampleMap 'primary'
splicing_CCP_P_plot <- left_join(splice_genes_mean, 
          colData(multiomics_MAE) %>% as_tibble() %>% select(name=patient_ID, PG ), 
          by=c("name")) %>%
  filter(!is.na(PG)) %>%
  ggplot(aes( PG, value ))+
  geom_boxplot(aes(fill=PG)) + geom_beeswarm() + 
  ggtitle("Spliceosome protein abundance") + 
  pp_sra +
  ylab("Mean abundance of spliceosome proteins")+
  #stat_compare_means(label = "p.signif", method = "t.test",
  #                   ref.group = ".all.", label.y = 0.2, hide.ns = TRUE)+
  guides(fill=guide_legend(title="PG")) +
  scale_fill_manual(values = colors_CCP)
splicing_CCP_P_plot +   theme(aspect.ratio=1) + 
    stat_compare_means(method = "anova", label.y = 0.22, hjust=0)

tmp_splice <- wideFormat(multiomics_MAE[splice_genes, ,"proteomics"], colDataCols = c("PG", "IGHV_mutated")) %>% as_tibble()
## Warning: 'experiments' dropped; see 'metadata'
## harmonizing input:
##   removing 522 sampleMap rows not in names(experiments)
##   removing 13 colData rownames not in sampleMap 'primary'
tmp_splice_mx <- tmp_splice %>% dplyr::select(-primary, -PG, -IGHV_mutated) %>% as.matrix()
rownames(tmp_splice_mx) <- tmp_splice$primary
colnames(tmp_splice_mx) <- gsub("proteomics_", "", colnames(tmp_splice_mx))
tmp_splice_anno <- tmp_splice %>% dplyr::select(PG, IGHV_mutated) %>% as.data.frame()
rownames(tmp_splice_anno) <- tmp_splice$primary
tmp_splice_anno$PG <- factor(tmp_splice_anno$PG,levels = c(5,1:4,6))

message("Heatmap of abundance spliceosome proteins ordered according to PG")
## Heatmap of abundance spliceosome proteins ordered according to PG
pheatmap(tmp_splice_mx[(tmp_splice_anno %>% rownames_to_column() %>% arrange( PG, IGHV_mutated) %>% .$rowname), ], 
         annotation_row = tmp_splice_anno, annotation_colors = ann_colors, scale = "column", cluster_rows = FALSE, color = inferno(10), border_color = NA,
         gaps_row = (which(!tmp_splice_anno %>% rownames_to_column() %>% arrange(PG, IGHV_mutated) %>% .$PG %>% duplicated())-1)[-1])

Selected spliceosome proteins

some_spliceprots <- c("SF3B1", "SNRPA", "PRPF6", "PRPF3", "SF3A1", "SNRPD2",  "SRSF4", "CDC5L", "PRPF19",
                      "CRNKL1", "PUF60", "PRPF8",
                      "SNRPB2")

tmp_Splice_GO <- wideFormat(multiomics_MAE[some_spliceprots, ,"proteomics"], colDataCols = c("PG", "IGHV_mutated" )) %>% 
  as_tibble()
## Warning: 'experiments' dropped; see 'metadata'
## harmonizing input:
##   removing 522 sampleMap rows not in names(experiments)
##   removing 13 colData rownames not in sampleMap 'primary'
tmp_Splice_GO_mx <- tmp_Splice_GO %>% dplyr::select(-primary, -PG, -IGHV_mutated) %>% as.matrix()
rownames(tmp_Splice_GO_mx) <- tmp_Splice_GO$primary
colnames(tmp_Splice_GO_mx) <- gsub("proteomics_", "", colnames(tmp_Splice_GO_mx))
tmp_Splice_GO_anno <- tmp_Splice_GO[, c("PG", "IGHV_mutated")] %>% as.data.frame()
rownames(tmp_Splice_GO_anno) <- tmp_Splice_GO$primary
tmp_Splice_GO_anno$PG <- factor(tmp_Splice_GO_anno$PG,levels = c(5,1:4,6))
SF3B1mut <- wideFormat(multiomics_MAE["SF3B1",,"SNPs"])$SNPs_SF3B1
## Warning: 'experiments' dropped; see 'metadata'
## harmonizing input:
##   removing 510 sampleMap rows not in names(experiments)
##   removing 1 colData rownames not in sampleMap 'primary'
names(SF3B1mut) <- wideFormat(multiomics_MAE["SF3B1",,"SNPs"])$primary
## Warning: 'experiments' dropped; see 'metadata'
## harmonizing input:
##   removing 510 sampleMap rows not in names(experiments)
##   removing 1 colData rownames not in sampleMap 'primary'
tmp_Splice_GO_anno$SF3B1 <- SF3B1mut[rownames(tmp_Splice_GO_anno)]
tmp_Splice_GO_anno$SF3B1[tmp_Splice_GO_anno$SF3B1 == 1] <- "mut"
tmp_Splice_GO_anno$SF3B1[tmp_Splice_GO_anno$SF3B1==0] <- "wt"
tmp_Splice_GO_anno$SF3B1 <- as.factor(tmp_Splice_GO_anno$SF3B1)
tmp_Splice_GO_anno$IGHV[tmp_Splice_GO_anno$IGHV_mutated==0] <- "U-CLL"
tmp_Splice_GO_anno$IGHV[tmp_Splice_GO_anno$IGHV_mutated==1] <- "M-CLL"
tmp_Splice_GO_anno <- tmp_Splice_GO_anno %>% select(-IGHV_mutated)

ann_colors = list(
  PG=c("1"= colors_CCP[1], "2"= colors_CCP[2], "3"= colors_CCP[3], "4"= colors_CCP[4], "5"= colors_CCP[5], "6"= colors_CCP[6] ),
  SF3B1=c("mut"="darkblue", "wt"="gray80"),
  IGHV=c("U-CLL"= "#0571b0", "M-CLL"= "#ca0020"))

pat_order_hclust_splice <- sapply(c(5,1:4,6), function(P){
  hc <- hclust(dist(tmp_Splice_GO_mx[rownames(tmp_Splice_GO_anno[tmp_Splice_GO_anno$PG==P,]), ] ))
  hc$labels[hc$order]
}) %>% unlist

PG_splice_proteins_pheatmap <- pheatmap(t(tmp_Splice_GO_mx[pat_order_hclust_splice, ]), 
         annotation_col = tmp_Splice_GO_anno, annotation_colors = ann_colors, scale = "row", cluster_cols = FALSE, 
         color =  inferno(length(breaks)), border_color = NA, 
         gaps_col = (which(!tmp_Splice_GO_anno %>% rownames_to_column() %>% arrange(PG, IGHV) %>% .$PG %>% duplicated())-1)[-1],
         breaks = breaks , show_colnames = F , treeheight_row = 0, fontsize_row = 8)

PG_splice_proteins_pheatmap_notscaled <- pheatmap(t(tmp_Splice_GO_mx[pat_order_hclust_splice, ]), 
         annotation_col = tmp_Splice_GO_anno, annotation_colors = ann_colors, cluster_cols = FALSE, 
         color =  inferno(length(breaks)), border_color = NA, 
         gaps_col = (which(!tmp_Splice_GO_anno %>% rownames_to_column() %>% arrange(PG, IGHV) %>% .$PG %>% duplicated())-1)[-1],
          show_colnames = F , treeheight_row = 0, fontsize_row = 8)

Spliceosome abundance across groups with SF3B1 mutational status color coded

left_join(splice_genes_mean, 
          wideFormat(multiomics_MAE["SF3B1", , "SNPs"], colDataCols = "PG") %>% as_tibble(), 
          by=c("name"="primary")) %>%
  filter(!is.na(PG)) %>%
  ggplot(aes( PG, value ))+
  geom_boxplot(aes(fill=PG)) + geom_beeswarm(aes(color=as.factor(SNPs_SF3B1))) + 
  stat_compare_means(method = "anova", label.y = 0.22, hjust=0) + 
  ggtitle("Spliceosome protein abundance") + 
  pp_sra_nox +
  ylab("Mean abundance of spliceosome proteins")+
  theme(aspect.ratio=1) +
  stat_compare_means(label = "p.signif", method = "t.test",
                     ref.group = ".all.", label.y = 0.2, hide.ns = TRUE)+
  guides(fill=guide_legend(title="PG"))+
  scale_fill_manual(values = colors_CCP) +
  scale_color_manual(values = c("grey", "red"))
## Warning: 'experiments' dropped; see 'metadata'
## harmonizing input:
##   removing 510 sampleMap rows not in names(experiments)
##   removing 1 colData rownames not in sampleMap 'primary'
## Warning: `gather_()` was deprecated in tidyr 1.2.0.
## Please use `gather()` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was generated.

BCAA

KEGG Valine, leucine and isoleucine degradation

BCAA_genes_mean <- assay(multiomics_MAE[BCAA_genes, ,"proteomics"]) %>% colMeans(na.rm = TRUE) %>% enframe()
## Warning: 'experiments' dropped; see 'metadata'
## harmonizing input:
##   removing 522 sampleMap rows not in names(experiments)
##   removing 13 colData rownames not in sampleMap 'primary'
BCAA_CCP_P_plot <- left_join(BCAA_genes_mean, 
          colData(multiomics_MAE) %>% as_tibble() %>% select(name=patient_ID, PG ), 
          by=c("name")) %>%
  filter(!is.na(PG)) %>%
  ggplot(aes( PG, value ))+
  geom_boxplot(aes(fill=PG)) + geom_beeswarm() + 
  ggtitle("BCAA protein abundance") + 
  pp_sra +
  ylab("Mean abundance of BCAA proteins")+
  #stat_compare_means(label = "p.signif", method = "t.test",
  #                   ref.group = ".all.", label.y = 0.55, hide.ns = TRUE)+
  guides(fill=guide_legend(title="PG")) +
  scale_fill_manual(values = colors_CCP)
BCAA_CCP_P_plot +   theme(aspect.ratio=1) + 
    stat_compare_means(method = "anova", label.y = 0.57, hjust=0)

Heatmap Selected proteins

some_BCAA <- c("ACAT1", "PCCA", "PCCB", "HADH", "HADHB", "HADHA", "HIBADH", "MCCC1", "MCCC2", "MLYCD")

tmp_BCAA <- wideFormat(multiomics_MAE[some_BCAA, ,"proteomics"], colDataCols = c("PG", "IGHV_mutated" )) %>% 
  as_tibble()
## Warning: 'experiments' dropped; see 'metadata'
## harmonizing input:
##   removing 522 sampleMap rows not in names(experiments)
##   removing 13 colData rownames not in sampleMap 'primary'
tmp_BCAA_mx <- tmp_BCAA %>% dplyr::select(-primary, -PG, -IGHV_mutated) %>% as.matrix()
rownames(tmp_BCAA_mx) <- tmp_BCAA$primary
colnames(tmp_BCAA_mx) <- gsub("proteomics_", "", colnames(tmp_BCAA_mx))
tmp_BCAA_anno <- tmp_BCAA[, c("PG", "IGHV_mutated")] %>% as.data.frame()
rownames(tmp_BCAA_anno) <- tmp_BCAA$primary
tmp_BCAA_anno$PG <- factor(tmp_BCAA_anno$PG,levels = c(5,1:4,6))
tmp_BCAA_anno$IGHV[tmp_BCAA_anno$IGHV_mutated==0] <- "U-CLL"
tmp_BCAA_anno$IGHV[tmp_BCAA_anno$IGHV_mutated==1] <- "M-CLL"
tmp_BCAA_anno <- tmp_BCAA_anno %>% select(-IGHV_mutated)

breaks= seq(min(tmp_BCAA_mx), max(tmp_BCAA_mx), 0.1)^2
breaks= sort(c(-breaks, breaks))
breaks <- breaks[! (breaks < min(tmp_BCAA_mx) | breaks > max(tmp_BCAA_mx) )]

ann_colors = list(
  PG=c("1"= colors_CCP[1], "2"= colors_CCP[2], "3"= colors_CCP[3], "4"= colors_CCP[4], "5"= colors_CCP[5], "6"= colors_CCP[6] ),
  IGHV=c("U-CLL"= "#0571b0", "M-CLL"= "#ca0020"))

pat_order_hclust_splice <- sapply(c(5,1:4,6), function(P){
  hc <- hclust(dist(tmp_BCAA_mx[rownames(tmp_BCAA_anno[tmp_BCAA_anno$PG==P,]), ] ))
  hc$labels[hc$order]
}) %>% unlist

PG_BCAA_proteins_pheatmap <- pheatmap(t(tmp_BCAA_mx[pat_order_hclust_splice, ]), 
         annotation_col = tmp_BCAA_anno, annotation_colors = ann_colors, scale = "row", cluster_cols = FALSE, 
         color =  inferno(length(breaks)), border_color = NA, 
         gaps_col = (which(!tmp_BCAA_anno %>% rownames_to_column() %>% arrange(PG, IGHV) %>% .$PG %>% duplicated())-1)[-1],
         breaks = breaks , show_colnames = F , treeheight_row = 0, fontsize_row = 5)

Proteasome

KEGG Proteoasome

proteasome_genes_mean <- assay(multiomics_MAE[proteasome_genes, ,"proteomics"]) %>% colMeans(na.rm = TRUE) %>% enframe()
## Warning: 'experiments' dropped; see 'metadata'
## harmonizing input:
##   removing 522 sampleMap rows not in names(experiments)
##   removing 13 colData rownames not in sampleMap 'primary'
proteasome_CCP_P_plot <- left_join(proteasome_genes_mean, 
          colData(multiomics_MAE) %>% as_tibble() %>% select(name=patient_ID, PG ), 
          by=c("name")) %>%
  filter(!is.na(PG)) %>%
  ggplot(aes( PG, value ))+
  geom_boxplot(aes(fill=PG)) + geom_beeswarm() + 
  ggtitle("Proteasome protein abundance") + 
  pp_sra +
  ylab("Mean abundance of proteasomal proteins")+
  #stat_compare_means(label = "p.signif", method = "t.test",
  #                   ref.group = ".all.", label.y = 0.3, hide.ns = TRUE)+
  guides(fill=guide_legend(title="PG")) +
  scale_fill_manual(values = colors_CCP)
proteasome_CCP_P_plot +   theme(aspect.ratio=1) + 
    stat_compare_means(method = "anova", label.y = 0.32, hjust=0)

Heatmap all proteasome proteins

some_proteasome <- proteasome_genes

tmp_proteasome <- wideFormat(multiomics_MAE[some_proteasome, ,"proteomics"], colDataCols = c("PG", "IGHV_mutated" )) %>% 
  as_tibble()
## Warning: 'experiments' dropped; see 'metadata'
## harmonizing input:
##   removing 522 sampleMap rows not in names(experiments)
##   removing 13 colData rownames not in sampleMap 'primary'
tmp_proteasome_mx <- tmp_proteasome %>% dplyr::select(-primary, -PG, -IGHV_mutated) %>% as.matrix()
rownames(tmp_proteasome_mx) <- tmp_proteasome$primary
colnames(tmp_proteasome_mx) <- gsub("proteomics_", "", colnames(tmp_proteasome_mx))
tmp_proteasome_anno <- tmp_proteasome[, c("PG", "IGHV_mutated")] %>% as.data.frame()
rownames(tmp_proteasome_anno) <- tmp_proteasome$primary
tmp_proteasome_anno$PG <- factor(tmp_proteasome_anno$PG,levels = c(5,1:4,6))
tmp_proteasome_anno$IGHV[tmp_proteasome_anno$IGHV_mutated==0] <- "U-CLL"
tmp_proteasome_anno$IGHV[tmp_proteasome_anno$IGHV_mutated==1] <- "M-CLL"
tmp_proteasome_anno <- tmp_proteasome_anno %>% select(-IGHV_mutated)

breaks= seq(min(tmp_proteasome_mx), max(tmp_proteasome_mx), 0.1)^2
breaks= sort(c(-breaks, breaks))
breaks <- breaks[! (breaks < min(tmp_proteasome_mx) | breaks > max(tmp_proteasome_mx) )]

ann_colors = list(
  PG=c("1"= colors_CCP[1], "2"= colors_CCP[2], "3"= colors_CCP[3], "4"= colors_CCP[4], "5"= colors_CCP[5], "6"= colors_CCP[6] ),
  IGHV=c("U-CLL"= "#0571b0", "M-CLL"= "#ca0020"))

pat_order_hclust_splice <- sapply(c(5,1:4,6), function(P){
  hc <- hclust(dist(tmp_proteasome_mx[rownames(tmp_proteasome_anno[tmp_proteasome_anno$PG==P,]), ] ))
  hc$labels[hc$order]
}) %>% unlist

PG_proteasome_proteins_pheatmap <- pheatmap(t(tmp_proteasome_mx[pat_order_hclust_splice, ]), 
         annotation_col = tmp_proteasome_anno, annotation_colors = ann_colors, scale = "row", cluster_cols = FALSE, 
         color =  inferno(length(breaks)), border_color = NA, 
         gaps_col = (which(!tmp_proteasome_anno %>% rownames_to_column() %>% arrange(PG, IGHV) %>% .$PG %>% duplicated())-1)[-1],
         breaks = breaks , show_colnames = F , treeheight_row = 0, fontsize_row = 5)

PCA

prot_pca <- prcomp(t( assay(multiomics_MAE[prot_few_nas , ,"proteomics"]) ))
## Warning: 'experiments' dropped; see 'metadata'
## harmonizing input:
##   removing 522 sampleMap rows not in names(experiments)
##   removing 13 colData rownames not in sampleMap 'primary'
summary(prot_pca)
## Importance of components:
##                           PC1     PC2     PC3     PC4    PC5     PC6     PC7
## Standard deviation     7.8245 6.66337 6.05805 5.39307 4.5962 4.47779 4.28476
## Proportion of Variance 0.1066 0.07734 0.06392 0.05066 0.0368 0.03492 0.03198
## Cumulative Proportion  0.1066 0.18398 0.24790 0.29856 0.3354 0.37028 0.40226
##                            PC8     PC9    PC10    PC11    PC12    PC13    PC14
## Standard deviation     4.07217 3.96780 3.89012 3.60729 3.50148 3.39828 3.33271
## Proportion of Variance 0.02888 0.02742 0.02636 0.02267 0.02136 0.02011 0.01935
## Cumulative Proportion  0.43114 0.45856 0.48492 0.50759 0.52894 0.54906 0.56840
##                           PC15    PC16    PC17   PC18    PC19   PC20    PC21
## Standard deviation     3.26312 3.19537 3.11575 2.9922 2.91051 2.8649 2.82250
## Proportion of Variance 0.01855 0.01778 0.01691 0.0156 0.01475 0.0143 0.01388
## Cumulative Proportion  0.58695 0.60473 0.62164 0.6372 0.65199 0.6663 0.68017
##                          PC22    PC23    PC24    PC25    PC26    PC27   PC28
## Standard deviation     2.7216 2.62868 2.61514 2.56747 2.53808 2.47949 2.4670
## Proportion of Variance 0.0129 0.01204 0.01191 0.01148 0.01122 0.01071 0.0106
## Cumulative Proportion  0.6931 0.70510 0.71702 0.72850 0.73972 0.75043 0.7610
##                           PC29    PC30    PC31    PC32    PC33    PC34    PC35
## Standard deviation     2.42571 2.36529 2.34684 2.29926 2.28396 2.22802 2.19890
## Proportion of Variance 0.01025 0.00974 0.00959 0.00921 0.00909 0.00865 0.00842
## Cumulative Proportion  0.77128 0.78102 0.79061 0.79982 0.80891 0.81755 0.82598
##                           PC36    PC37    PC38    PC39    PC40    PC41    PC42
## Standard deviation     2.15073 2.13611 2.08023 2.05930 2.03908 2.01421 1.99135
## Proportion of Variance 0.00806 0.00795 0.00754 0.00739 0.00724 0.00707 0.00691
## Cumulative Proportion  0.83403 0.84198 0.84952 0.85690 0.86415 0.87121 0.87812
##                           PC43   PC44    PC45    PC46    PC47    PC48    PC49
## Standard deviation     1.97767 1.9325 1.90702 1.86844 1.84392 1.83474 1.81373
## Proportion of Variance 0.00681 0.0065 0.00633 0.00608 0.00592 0.00586 0.00573
## Cumulative Proportion  0.88493 0.8914 0.89777 0.90385 0.90978 0.91564 0.92137
##                           PC50    PC51    PC52    PC53    PC54   PC55    PC56
## Standard deviation     1.77205 1.76322 1.74185 1.72077 1.70172 1.6595 1.63090
## Proportion of Variance 0.00547 0.00542 0.00528 0.00516 0.00504 0.0048 0.00463
## Cumulative Proportion  0.92684 0.93225 0.93754 0.94270 0.94774 0.9525 0.95717
##                           PC57    PC58    PC59    PC60   PC61    PC62    PC63
## Standard deviation     1.61852 1.59146 1.56141 1.53137 1.5157 1.49844 1.46325
## Proportion of Variance 0.00456 0.00441 0.00425 0.00408 0.0040 0.00391 0.00373
## Cumulative Proportion  0.96173 0.96614 0.97039 0.97447 0.9785 0.98239 0.98612
##                           PC64    PC65    PC66   PC67      PC68
## Standard deviation     1.44810 1.43129 1.39002 1.3758 6.207e-15
## Proportion of Variance 0.00365 0.00357 0.00337 0.0033 0.000e+00
## Cumulative Proportion  0.98977 0.99334 0.99670 1.0000 1.000e+00
prot_pca_x <- as_tibble(prot_pca$x[,1:10])
prot_pca_x$pat_ID <- colnames(assay(multiomics_MAE[prot_few_nas , ,"proteomics"]))
## Warning: 'experiments' dropped; see 'metadata'
## harmonizing input:
##   removing 522 sampleMap rows not in names(experiments)
##   removing 13 colData rownames not in sampleMap 'primary'
prot_pca_x <- left_join(prot_pca_x, 
                        wideFormat(multiomics_MAE[, ,c("SNPs","chrom_abber", "health_record_bin") ], 
                                   colDataCols = c("gender", "treatment_status", "PG")) %>% as_tibble() , 
                        by=c("pat_ID"="primary") )
## Warning: 'experiments' dropped; see 'metadata'
## Warning: 'ExperimentList' contains 'data.frame' or 'DataFrame',
##   potential for errors with mixed data types
## harmonizing input:
##   removing 348 sampleMap rows not in names(experiments)
#prot_pca_x <- prot_pca_x %>% replace(is.na(.), "unknown")

ggplot(prot_pca_x, aes(PC1, PC2)) +
  geom_point(aes(color=as.factor(health_record_bin_IGHV_mutated), shape=as.factor(chrom_abber_trisomy12)))  + pp_sra +
  scale_shape_manual(values = c( 16,  1, 4)) + 
  scale_color_manual(values = c( "#0571b0","#ca0020", "grey")) + theme(aspect.ratio = 1)
## Warning: Removed 4 rows containing missing values (geom_point).

ggplot(prot_pca_x, aes(PC1, PC3)) +geom_point(aes(color=as.factor(health_record_bin_IGHV_mutated), shape=as.factor(chrom_abber_trisomy12)) ) +
  scale_shape_manual(values = c( 16,  1, 4)) + 
  scale_color_manual(values = c( "#0571b0","#ca0020", "grey"))+ pp_sra+ theme(aspect.ratio = 1)
## Warning: Removed 4 rows containing missing values (geom_point).

ggplot(prot_pca_x, aes(PC2, PC3)) +geom_point(aes(color=as.factor(health_record_bin_IGHV_mutated), shape=as.factor(chrom_abber_trisomy12)) ) +  
  scale_shape_manual(values = c( 16,  1, 4)) + 
  scale_color_manual(values = c( "#0571b0", "#ca0020", "grey")) + pp_sra+ theme(aspect.ratio = 1)
## Warning: Removed 4 rows containing missing values (geom_point).

ggplot(prot_pca_x, aes(PC1, PC2)) +geom_point(aes(color=gender)) + pp_sra + scale_color_manual(values = c( "#0571b0","#ca0020", "grey")) + theme(aspect.ratio = 1)

PCA_CCP_1_2 <- ggplot(prot_pca_x, aes(PC1, PC2)) +
  geom_point(aes(fill=PG),shape=21, color="grey") + pp_sra  +
  guides(color=guide_legend(title="PG")) +
  scale_fill_manual(values = colors_CCP) +
  xlab( paste0("PC1 (", round(summary(prot_pca)$importance[2,1], 2)*100, "%)" ))+
  ylab( paste0("PC2 (", round(summary(prot_pca)$importance[2,2], 2)*100, "%)" ))
PCA_CCP_1_2 +
  theme(aspect.ratio=1, plot.title = element_text(size = 30))

PCA_CCP_1_3 <- ggplot(prot_pca_x, aes(PC1, PC3)) +
  geom_point(aes(fill=PG),shape=21, color="grey") + pp_sra  +
  guides(color=guide_legend(title="PG")) +
  scale_fill_manual(values = colors_CCP)  +
  xlab( paste0("PC1 (", round(summary(prot_pca)$importance[2,1], 2)*100, "%)" ))+
  ylab( paste0("PC3 (", round(summary(prot_pca)$importance[2,3], 2)*100, "%)" ))
PCA_CCP_1_3 +
  theme(aspect.ratio=1, plot.title = element_text(size = 30))

PCA_CCP_2_3 <- ggplot(prot_pca_x, aes(PC2, PC3)) +
  geom_point(aes(fill=PG),shape=21, color="grey") + pp_sra  +
  guides(color=guide_legend(title="PG")) +
  scale_fill_manual(values = colors_CCP)  +
  xlab( paste0("PC2 (", round(summary(prot_pca)$importance[2,2], 2)*100, "%)" ))+
  ylab( paste0("PC3 (", round(summary(prot_pca)$importance[2,3], 2)*100, "%)" ))
PCA_CCP_2_3 +
  theme(aspect.ratio=1, plot.title = element_text(size = 30))

PCA_CCP_1_4 <- ggplot(prot_pca_x, aes(PC1, PC4)) +
  geom_point(aes(fill=PG),shape=21, color="grey") + pp_sra  +
  guides(color=guide_legend(title="PG")) +
  scale_fill_manual(values = colors_CCP)  +
  xlab( paste0("PC1 (", round(summary(prot_pca)$importance[2,1], 2)*100, "%)" ))+
  ylab( paste0("PC4 (", round(summary(prot_pca)$importance[2,4], 2)*100, "%)" ))
PCA_CCP_1_4 +
  theme(aspect.ratio=1, plot.title = element_text(size = 30))

tsne RNA

genes_no_nas <- multiomics_MAE[["RNAseq_norm"]] %>% is.na() %>% rowSums()
genes_no_nas <- genes_no_nas[ genes_no_nas == 0 ] %>% names()

set.seed(10)
rtsne_out_RNA <- Rtsne( t(assay(multiomics_MAE[genes_no_nas , ,"RNAseq_norm"])), perplexity = 10 )
## Warning: 'experiments' dropped; see 'metadata'
## Warning: 'ExperimentList' contains 'data.frame' or 'DataFrame',
##   potential for errors with mixed data types

## Warning: 'ExperimentList' contains 'data.frame' or 'DataFrame',
##   potential for errors with mixed data types
rtsne_out_RNA_df <- rtsne_out_RNA$Y %>% as.data.frame() %>% as_tibble()
rtsne_out_RNA_df$pat_ID <- colnames(assay(multiomics_MAE[genes_no_nas , ,"RNAseq_norm"]))
## Warning: 'experiments' dropped; see 'metadata'

## Warning: 'ExperimentList' contains 'data.frame' or 'DataFrame',
##   potential for errors with mixed data types

## Warning: 'ExperimentList' contains 'data.frame' or 'DataFrame',
##   potential for errors with mixed data types
rtsne_out_RNA_df <- left_join(rtsne_out_RNA_df, 
                          (wideFormat(multiomics_MAE[, ,c("SNPs","chrom_abber", "health_record_bin") ], colDataCols = c("gender", "treatment_status", "doehner_groups", "PG", "CCP6_RNA")) %>% as_tibble()), 
                          by=c("pat_ID"="primary") )
## Warning: 'experiments' dropped; see 'metadata'

## Warning: 'ExperimentList' contains 'data.frame' or 'DataFrame',
##   potential for errors with mixed data types
rtsne_out_RNA_df <- mutate_at(rtsne_out_RNA_df, colnames(rtsne_out_RNA_df %>% dplyr::select(SNPs_ATM:health_record_bin_treated)), as.logical)

#rtsne_out_RNA_df <- rtsne_out_RNA_df %>% replace(is.na(.), "unknown")

tsne_CCP_R_plot <- ggplot(rtsne_out_RNA_df, aes(V1, V2)) +
  geom_point(aes(fill=PG), color="grey", shape=21) +
  scale_fill_manual(values = colors_CCP)+
  pp_sra +
  guides(color=guide_legend(title="PG"))
tsne_CCP_R_plot + theme(aspect.ratio=1, plot.title = element_text(size = 30))

tsne_CCPRNA_R_plot <- ggplot(rtsne_out_RNA_df, aes(V1, V2)) +
  geom_point(aes(fill=CCP6_RNA), color="grey", shape=21) +
  scale_color_hue()+
  pp_sra +
  guides(color=guide_legend(title="CC_RNA"))
tsne_CCPRNA_R_plot+  theme(aspect.ratio=1, plot.title = element_text(size = 30))

sapply(colnames(rtsne_out_RNA_df)[c(4,5,9:16)], function(var){
  print(ggplot(rtsne_out_RNA_df, aes(V1, V2)) +
          geom_point(aes_string(color=var), size=3) +
    scale_color_manual(values = c("#92c5de", "#f4a582",  "grey"))  + 
          ggtitle(gsub("SNPs_", "", gsub("chrom_abber_", "", gsub( "health_record_bin_", "", var) ))) +
      pp_sra_noguides +
      theme(aspect.ratio=1, plot.title = element_text(size = 30)) )
})

##             gender           treatment_status SNPs_ATM         SNPs_BIRC3      
## data        tbl_df,27        tbl_df,27        tbl_df,27        tbl_df,27       
## layers      list,1           list,1           list,1           list,1          
## scales      ScalesList,2     ScalesList,2     ScalesList,2     ScalesList,2    
## mapping     uneval,2         uneval,2         uneval,2         uneval,2        
## theme       theme,93         theme,93         theme,93         theme,93        
## coordinates CoordCartesian,5 CoordCartesian,5 CoordCartesian,5 CoordCartesian,5
## facet       FacetNull,2      FacetNull,2      FacetNull,2      FacetNull,2     
## plot_env    ?                ?                ?                ?               
## labels      list,4           list,4           list,4           list,4          
## guides      list,4           list,4           list,4           list,4          
##             SNPs_EGR2        SNPs_NOTCH1      SNPs_POT1        SNPs_SF3B1      
## data        tbl_df,27        tbl_df,27        tbl_df,27        tbl_df,27       
## layers      list,1           list,1           list,1           list,1          
## scales      ScalesList,2     ScalesList,2     ScalesList,2     ScalesList,2    
## mapping     uneval,2         uneval,2         uneval,2         uneval,2        
## theme       theme,93         theme,93         theme,93         theme,93        
## coordinates CoordCartesian,5 CoordCartesian,5 CoordCartesian,5 CoordCartesian,5
## facet       FacetNull,2      FacetNull,2      FacetNull,2      FacetNull,2     
## plot_env    ?                ?                ?                ?               
## labels      list,4           list,4           list,4           list,4          
## guides      list,4           list,4           list,4           list,4          
##             SNPs_TP53        SNPs_XPO1       
## data        tbl_df,27        tbl_df,27       
## layers      list,1           list,1          
## scales      ScalesList,2     ScalesList,2    
## mapping     uneval,2         uneval,2        
## theme       theme,93         theme,93        
## coordinates CoordCartesian,5 CoordCartesian,5
## facet       FacetNull,2      FacetNull,2     
## plot_env    ?                ?               
## labels      list,4           list,4          
## guides      list,4           list,4

Color tsne based on expression of B-cell receptor signaling pathways

BCR_genes_RNA <- left_join(enframe(BCR_genes, value = "symbol", name = NULL), metadata(multiomics_MAE)$gene_symbol_mapping, by=c("symbol"="hgnc_symbol")) %>% 
  filter(!is.na(ensembl_gene_id)) %>% .$ensembl_gene_id
BCR_genes_RNA_mean <- assay(multiomics_MAE[BCR_genes_RNA, ,"RNAseq_norm"]) %>% colMeans(na.rm = TRUE) %>% enframe()
## Warning: 'experiments' dropped; see 'metadata'
## Warning: 'ExperimentList' contains 'data.frame' or 'DataFrame',
##   potential for errors with mixed data types

## Warning: 'ExperimentList' contains 'data.frame' or 'DataFrame',
##   potential for errors with mixed data types

Save important plots

save(tsne_CCP_P_plot, tsne_CCP_R_plot,tsne_CCPRNA_R_plot, tsne_IGHVtris_plot, BCR_CCP_P_plot, splicing_CCP_P_plot,
     PCA_CCP_1_2, PCA_CCP_1_3, PCA_CCP_1_4, PCA_CCP_2_3, 
     PG_BCR_proteins_pheatmap, PG_BCR_proteins_pheatmap_notscaled,
     PG_splice_proteins_pheatmap, PG_splice_proteins_pheatmap_notscaled, PG_BCAA_proteins_pheatmap,
     BCAA_CCP_P_plot, proteasome_CCP_P_plot, PG_proteasome_proteins_pheatmap,
file = "RData_plots/CLL_Proteomics_DimensionReduction_Plots.RData")

Session Info

sessionInfo()
## R version 4.1.3 (2022-03-10)
## Platform: x86_64-apple-darwin17.0 (64-bit)
## Running under: macOS Big Sur/Monterey 10.16
## 
## Matrix products: default
## BLAS:   /Library/Frameworks/R.framework/Versions/4.1/Resources/lib/libRblas.0.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.1/Resources/lib/libRlapack.dylib
## 
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
## 
## attached base packages:
## [1] stats4    stats     graphics  grDevices utils     datasets  methods  
## [8] base     
## 
## other attached packages:
##  [1] viridis_0.6.2               viridisLite_0.4.0          
##  [3] gplots_3.1.3                MultiAssayExperiment_1.20.0
##  [5] Rtsne_0.16                  limma_3.50.3               
##  [7] fdrtool_1.2.17              vsn_3.62.0                 
##  [9] forcats_0.5.1               stringr_1.4.0              
## [11] dplyr_1.0.9                 purrr_0.3.4                
## [13] readr_2.1.2                 tidyr_1.2.0                
## [15] tibble_3.1.7                tidyverse_1.3.1            
## [17] DESeq2_1.34.0               SummarizedExperiment_1.24.0
## [19] Biobase_2.54.0              MatrixGenerics_1.6.0       
## [21] matrixStats_0.62.0          GenomicRanges_1.46.1       
## [23] GenomeInfoDb_1.30.1         IRanges_2.28.0             
## [25] S4Vectors_0.32.4            BiocGenerics_0.40.0        
## [27] ggbeeswarm_0.6.0            ggpubr_0.4.0               
## [29] Hmisc_4.7-0                 ggplot2_3.3.6              
## [31] Formula_1.2-4               survival_3.3-1             
## [33] lattice_0.20-45             Matrix_1.4-1               
## [35] pheatmap_1.0.12             gtools_3.9.3               
## [37] plyr_1.8.7                 
## 
## loaded via a namespace (and not attached):
##   [1] readxl_1.4.0           backports_1.4.1        splines_4.1.3         
##   [4] BiocParallel_1.28.3    digest_0.6.29          htmltools_0.5.2       
##   [7] fansi_1.0.3            magrittr_2.0.3         checkmate_2.1.0       
##  [10] memoise_2.0.1          cluster_2.1.3          tzdb_0.3.0            
##  [13] Biostrings_2.62.0      annotate_1.72.0        modelr_0.1.8          
##  [16] jpeg_0.1-9             colorspace_2.0-3       blob_1.2.3            
##  [19] rvest_1.0.2            haven_2.5.0            xfun_0.31             
##  [22] crayon_1.5.1           RCurl_1.98-1.7         jsonlite_1.8.0        
##  [25] genefilter_1.76.0      glue_1.6.2             gtable_0.3.0          
##  [28] zlibbioc_1.40.0        XVector_0.34.0         DelayedArray_0.20.0   
##  [31] car_3.1-0              abind_1.4-5            scales_1.2.0          
##  [34] DBI_1.1.3              rstatix_0.7.0          Rcpp_1.0.9            
##  [37] xtable_1.8-4           htmlTable_2.4.1        foreign_0.8-82        
##  [40] bit_4.0.4              preprocessCore_1.56.0  htmlwidgets_1.5.4     
##  [43] httr_1.4.3             RColorBrewer_1.1-3     ellipsis_0.3.2        
##  [46] farver_2.1.1           pkgconfig_2.0.3        XML_3.99-0.10         
##  [49] nnet_7.3-17            sass_0.4.1             dbplyr_2.2.1          
##  [52] deldir_1.0-6           locfit_1.5-9.6         utf8_1.2.2            
##  [55] labeling_0.4.2         tidyselect_1.1.2       rlang_1.0.4           
##  [58] AnnotationDbi_1.56.2   munsell_0.5.0          cellranger_1.1.0      
##  [61] tools_4.1.3            cachem_1.0.6           cli_3.3.0             
##  [64] generics_0.1.3         RSQLite_2.2.14         broom_1.0.0           
##  [67] evaluate_0.15          fastmap_1.1.0          yaml_2.3.5            
##  [70] knitr_1.39             bit64_4.0.5            fs_1.5.2              
##  [73] caTools_1.18.2         KEGGREST_1.34.0        xml2_1.3.3            
##  [76] compiler_4.1.3         rstudioapi_0.13        beeswarm_0.4.0        
##  [79] png_0.1-7              affyio_1.64.0          ggsignif_0.6.3        
##  [82] reprex_2.0.1           geneplotter_1.72.0     bslib_0.3.1           
##  [85] stringi_1.7.8          highr_0.9              vctrs_0.4.1           
##  [88] pillar_1.7.0           lifecycle_1.0.1        BiocManager_1.30.18   
##  [91] jquerylib_0.1.4        data.table_1.14.2      bitops_1.0-7          
##  [94] R6_2.5.1               latticeExtra_0.6-30    affy_1.72.0           
##  [97] KernSmooth_2.23-20     gridExtra_2.3          vipor_0.4.5           
## [100] assertthat_0.2.1       withr_2.5.0            GenomeInfoDbData_1.2.7
## [103] parallel_4.1.3         hms_1.1.1              grid_4.1.3            
## [106] rpart_4.1.16           rmarkdown_2.14         carData_3.0-5         
## [109] lubridate_1.8.0        base64enc_0.1-3        interp_1.1-3